package scu.maqiang.possion;

import scu.maqiang.fes.BVPType;
import scu.maqiang.fes.FES2T31;
import scu.maqiang.mesh.Mesh2T3;
import scu.maqiang.mesh.ScalarFunc;
import scu.maqiang.mesh.T3Type;
import scu.maqiang.numeric.MVO;

/**
 * 轴对称Possion问题
 */
public class AxsyPossion2DT31 {

	public static void main(String[] args) {
		Mesh2T3 mesh = new Mesh2T3().square2D(400, 40, T3Type.Right);
		mesh.scale(1.0, 0.1);
		mesh.translate(1.0, 0.0);
		
		Possion2DT31 p = new Possion2DT31(mesh);
		p.DirichletProblem(BVPType.AXISYMMETRIC, new double[]{1.0}, 10.0, new int[] {2, 4},  0.0, null);
		double[] x = p.exportSolution();
		
		FES2T31 fs = new FES2T31(mesh);
		ScalarFunc rf = (xy, label, param) -> {
			return -2.50 * xy[0] * xy[0] + 7.5 * Math.log(xy[0]) / Math.log(2.0) + 2.50;
		};
		System.out.println("Result min: " + MVO.min(x) + "\tmax: " + MVO.max(x));
		double[] xExact = fs.valueFromFunc(rf, null);
		System.out.println("Exact min: " + MVO.min(xExact) + "\tmax: " + MVO.max(xExact));
		double[] error0 = MVO.add(xExact, -1.0, x);
		System.out.println("Error min: " + MVO.min(error0) + "\tmax: " + MVO.max(error0));
		double L2Norm = fs.computeL2Norm(xExact);
		System.out.println("L2 Error: " + fs.computeL2Norm(error0) / L2Norm);
	}
}